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S" 

f"^ ■ Surface aging phenomena are discussed for semi-infinite systems prepared in a fully disordered 

initial state and then quenched to or below the critical point. Besides solving exactly the semi-infinite 
£N| , Ising model in the limit of large dimensions, we also present results of an extensive numerical study of 

the nonequilibrium dynamical behavior of the two-dimensional semi-infinite Ising model undergoing 
coarsening. The studied models reveal a simple aging behavior where some of the nonequilibrium 
surface exponents take on values that differ from their bulk counterparts. For the two-dimensional 
. semi-infinite Ising model we find that the exponent 61, that describes the scaling behavior of the 

' surface autocorrelation, vanishes. These simulations also reveal the existence of strong finite-time 

corrections that to some extent mask the leading scaling behavior of the studied two-time quantities. 
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I. INTRODUCTION 

Intriguing phenomena are observed when bringing simple fcrromagncts out of equilibrium through a temperature 
quenchii 2 ^!^. Consider a ferromagnetic system prepared at high temperatures, i.e. in a disordered and uncorrelated 
initial state, that is suddenly quenched onto or below the critical point. If the quench is onto the critical point, 

■ critical dynamics sets in, yielding a dynamical correlation length that increases with time. For a quench inside the 
1 I ordered phase, the formation and growth of well-ordered domains are observed. Interestingly, these two physically very 

■ different cases have in common that the relevant length scale L (which in the first case is the dynamical correlation 
length, whereas in the second case it is given by the typical extent of the ordered domains) increases as a simple 

O . power-law of time: 

L(t)~tV«. (1) 

^ ' Assuming non-conserved dynamics (which is the only dynamics studied in this paper) one finds for the dynamical 

OO , exponent z the value 2 below the critical temperature T c , whereas at T c the dynamical exponent may take on values 

,— ' ■ slightly larger than 2. 

£NJ ' The power-law growth ([1]) is responsible for many nonequilibrium phenomena which arc usually summarized under 
the header of simple aging. Thus, it follows directly from {!]) that two-time quantities, like dynamical correlation and 

■ response functions, display dynamical scaling. Introducing the time- and space-dependent order parameter <j)(r;t), 
f* 1 **. \ the correlation function can be written as 



C(t,s;r-r') = {<j>{r;t)<f>{r';s)), (2) 



whereas the response function, which measures the response of the order parameter at site r at time t to an external 
$H ' field h(f'\s) acting on site r ' at time s, is defined by 



5h{f';s) 



(t > s). (3) 

h=0 



In writing these equations we assume spatial translation invariancc, as encountered in ferromagnetic bulk systems. 
The usually studied autocorrelation function C(t, s) and autoresponse function R(t, s) are obtained by setting r = r ' 
in Equations © and ([3]), respectively. In the dynamical scaling regime with t, s, t — s 3> T m i cro , where T m icro is 
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a microscopic time scale, these two-time quantities can be cast into a simple scaling form. For example, for the 
autocorrelation and the autoresponse functions we have 

C(t, a) = s- h f c (t/s) and R(t, s) = s" 1 ^ f R (t/ S), (4) 

where a and b are noncquilibrium exponents, whereas fc and Jr are scaling functions that for large arguments display 
a power-law decay, 

fc(y) ~ y~ Xo/z and f R (y)~y- XR/z , (5) 

with the autocorrelation 6 * 7 - and the autoresponse exponents^ Ac and A^. At the critical point, the nonequilibrium 
exponents a and b can be expressed by known critical exponents, yielding a = b = (d — 2 + r\)j z where d is the 
dimensionality of the system and 77 is the usual static critical exponent governing the power-law decay of the spatial 
correlations at equilibrium. In the ordered phase one has b = and a = 1/z for systems with exponentially decaying 
static correlation s] 1 ! 9 ' 10 . At the critical point the autocorrelation exponent is related to the so-called initial slip 
exponent^. In addition, autocorrelation and autoresponse exponents can be shown to be identical in systems with 
short-ranged initial correlations and purely relaxational dynamic s] 1 ' 12 . 

This briefly described simple aging scenario has been studied very intensively in bulk systems, but for systems 
bounded by surfaces the investigation of local aging processes close to surfaces is only at its very beginning 3 - ' 13 ' 14 . The 
emerging interest in surface aging phenomena can be related to the increasing importance of surface dominated small 
systems in research and technology. Indeed, nonequilibrium processes are deeply affected by the presence of surfaces 
which can result in changes in the physical behavior even at macroscopic distances from the surface. In principle, 
these surface properties can be studied by X-ray scattering at grazing incidence. 

Looking at critical systems bounded by surfaces, it is well known that already the static critical behavior at a 
surface is different from the bulk critical behavior, yielding new sets of static surface critical exponent s] 15 ' 16 ' 17 ' 18 . In 
fact, one even encounters different surface universality classes for a given bulk universality class, depending on the 
value of the surface couplings or on the existence of a surface field. Looking at the dynamics, it has to be noted that 
in the case of purely diffusive dynamics the dynamical exponent z has the same universal value close to the surface 
as inside the bul k 19 ' 20 ' 21 . One then expects a surface aging behavior similar to the bulk behavior, but with local 
nonequilibrium exponents and scaling functions of local two-time functions that differ from the corresponding bulk 
quantities. 

In order to discuss the expected aging phenomenology close to a critical surface in more detail, let us consider 
an idealized semi- infinite lattice in d dimensions where we write the position vector r as r — (x,y). Here x is a 
{d — l)-dimcnsional vector parallel to the surface, whereas y labels the layers perpendicular to the surface (with y = 1 
being the surface layer). With this we obtain the following generalizations for the correlation and response functions: 

C(t,s;y,y',x-x') 
R(t,s;y,y',x- x") 

where we assumed spatial translation invariance in the directions parallel to the surface. For y, y' — ► 00 we recover 
the bulk quantities, whereas y = y' = 1 yields the surface correlation and response functions. Of special interest arc 
the surface autocorrelation and autoresponse functions with x = x ' that we write as C\(t, s) := C(t,s; 1, 1,0 ) and 
Ri(t, s) := R(t, s; 1, 1, ). For these surface quantities, the simple scaling forms 

= s-^fcAt/s), fedt/s) ~ (V*r Aci/z (6) 
Rl (t,s) = s-^fRAt/s), f Rl (t/s) ~ (t/s)- x ^ z 

are expected 3 -^ when t, s and also the difference t — s are large compared to some microscopic timcscalc. The scaling 
functions /d(t/s) and fR 1 (t/s) should again display a simple power-law behavior for large values of t/s. General 
scaling arguments^ allow to express the surface exponents appearing in ^ through other known exponents: 

oi = h = (d - 2 + rj\\)/ z and A Cl = A^ = A c + - 77 (7) 

where 7711 is the static exponent that governs the decay of the correlations parallel to the surface. As for bulk 
system s] 22 ' 23 , surface autocorrelation and autoresponse functions can be combined to yield the surface fluctuation- 
dissipation ratioi 3 - 



(4>(x,y;t)4>(x ,y,s)) 
5 (<p{x,y;t)} 
5h(x',y';s) h=0 



(8) 
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with a universal limit value 

Xf° = lim ( lim Xi(t,sV) (9) 

that characterizes the different dynamical surface universality classes^. The scaling picture ([6]) and the relations 
between the various nonequilibrium exponents have been verified by one of us through a numerical study of the out- 
of-equilibrium dynamics of various critical semi-infinite Ising models^. In addition, the critical semi-infinite Gaussian 
model^ and the critical semi-infinite spherical model^ were also found to display this simple aging scenario. 

Whereas at least some knowledge has accumulated in recent years on the local aging behavior close to critical 
surfaces, almost nothing is known on surface aging processes taking place in coarsening systems. In RefJ^ we have 
looked at the out-of-equilibrium dynamical behavior of the semi-infinite spherical model. For this special model we 
have verified the existence of dynamical scaling and simple aging close to surfaces for quenches inside the ordered 
phase. Surprisingly, the nonequilibrium exponent b\, describing the scaling behavior of the surface autocorrelation, 
was found to take on the value b\ = 1 , different from the standard value b = of the corresponding exponent in bulk 
systems undergoing phase-ordering. This result calls for a thorough investigation of surface aging phenomena in other 
semi-infinite systems with phase-ordering dynamics. 

In this paper we continue our study of local aging processes in bounded fcrromagncts. On the one hand, we discuss 
the semi-infinite short-range Ising model in the limit of high dimensions that can be solved exactly. On the other 
hand, we present results of extensive Monte Carlo simulations of the standard two-dimensional semi-infinite Ising 
model prepared at high temperatures and then quenched inside the ordered phase. These numerical results yield new 
and interesting insights into the local processes taking place in coarsening systems close to surfaces. All the systems 
studied have in common that the dynamical exponent takes on the value z = 2. 

The paper is organized as follows. In the next Section we compute scaling functions and nonequilibrium exponents 
in the exactly solvable semi-infinite model. Our numerical results obtained from simulations of the two-dimensional 
semi-infinite Ising model undergoing phase-ordering are then presented in Section III. Finally, in Section IV we draw 
our conclusions and summarize our results. 



II. QUENCHING SEMI-INFINITE SYSTEMS FROM HIGH TEMPERATURES: EXACT RESULTS 



Exactly solvable models are often quite unrealistic and even artificial. One of the reasons for nevertheless studying 
this kind of models in physics is to obtain a guidance for the development of a future more sophisticated theoretical 
approach. With this in mind, we discuss in the following the nonequilibrium dynamical behavior of the exactly 
solvable short-range semi-infinite Ising model in the limit of high dimensions that is prepared in an uncorrelated 
initial state with vanishing magnetization and then quenched below or at the critical point. As in that limit the 
model is mean-field like, we expect the same critical exponents and the same scaling functions (up to some numerical 
prefactors) as those found in the semi-infinite Gaussian model^. 

The out-of-equilibrium behavior of the bulk Ising model with nearest neighbor ferromagnetic interactions has 
recently been studied in the limit of a large number d of space dimensions^. Here we generalize the calculations of 
Garriga et al. to the semi-infinite case. 

Using a semi- infinite hypcrcube with lattice constant 1, the Hamiltonian of our model can be written in the very 
general form 

H= ~Yd ^ <TS ' lCTS '' 1 ~ <^5Z H Vx,yV3\y - ^^2^2<Jx,y<Jx.y+l (10) 

(x,x ') y>2 (x,x ') y>i x 

where the sum over (x, x ') indicates a sum over all nearest neighbor pairs lying in the same layer. The spins can 
take on the values ±1, and an additional field term can be added if needed. In writing (|10p we take into account the 
layered structure of the lattice and distinguish between nearest neighbor pairs lying in a layer parallel to the surface 
and nearest neighbor pairs belonging to different layers. As usual when dealing with semi-infinite systems^, we have 
introduced a different coupling constant J s for interactions between nearest neighbor spins located both in the surface 
layer. We will however restrict ourselves in the following to the special case J s = Jb = 1. On the one hand this yields 
in the limit d — > oo the critical temperature T c = 1 (where we set ks = 1), on the other hand we then encounter at 
the critical temperature the so-called ordinary transition^ where the bulk alone is critical. 
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The main difference between the present case and the model considered in Ref^ is of course the absence of spatial 
translation invariance in the direction perpendicular to the surface. Due to this, the time-dependent local fields that 
the spins experience are now layer-dependent, leading to the expressions 



\ x'(x) J 

H v {t) = h% x i{t) + i- ag, y+1 (t) + <T 2ltf _i(t) + J2 ^'At) for v¥=i, (ii) 

V S'(g) J 

where the sum over x '(x) is the sum over the in-plane nearest neighbor lattice sites x ' of x. We added in these equa- 
tions an external field h^ x * (t) needed for the computation of the response function (for the computation of the correla- 
tion function, h e g*(t) is of course set to zero). Using heat- bath dynamics, these local fields hg y (t) (which depend on the 
dimension d) appear in the flip rates, as each spin will flip independently with the rate (1 — &g,y(t) tanh(/ij i2/ (f)/T))/2. 
It is important to note that in absence of an external magnetic field the magnetization remains at any time at is 
initial value zero everywhere in the sample. 



A. The correlation function 



In their paper— Garriga et al. derived very general equations of motion for the one- and the two-time correlation 
functions C and C that can also be used in our case by plugging in the layer-dependent local fields (fTTj) . Recalling 
that we still have invariance for spatial translations parallel to the surface, we can write the following equations: 

d t C(t;y,y',x-x') = -2C(t;y,y',x - x ') + (At 3 , y (t)Aa s ,, y ,(t)) 

+ (Aosf ltf (t)At*/ lW /(f)} (12) 
d t C(t,s;y,y',x-x') = -C(t, s; y,y',x - x ') + (Ats iV (t)Aag^yi(s)) (13) 

where we use the notations Atg tV (t) := t&nh(hs ty (t) /T) — (ta,nh(hs, y {t) /T)} and Aag_ y (t) = cr 3y (t) — {<Js,y{t)) for the 
deviations from the averages. 

In the limit of large d we can develop tanh(/i^ J/ (t)/T) in l/d (recall that no external field is acting on the spins and 
that the local fields hg y are layer-dependent) which then yields the following expressions for the equations of motion: 
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d t C{t;y,y',S) = -2C(t;y,y',x) + J(c(i;y+ l,y',x) + C(t;y- l,y',x) (14) 
+ C(t; y, y' + 1, x) + C(t; y, y'-l,x)+2 ^C(t; y, y' , z)) + b{t; y, y', x) 

z(x) 

d t C(t,s;y,y',x) = -C(t,s;y,y',x) + ^(c(t,s;y + l,y',x) 



with 



C(t, s;y-l, y', x) + ^ C(t, s; y, y', z)j . (15) 

z{3) 

C(t,s;0,y',x)=0 = C(t,s;y,0,x) . (16) 



In writing these equations we exploit the spatial translation invariance parallel to the surface by setting x ' = 0. The 
parameter 7 is given by 7 := l/(Td), whereas the sum over z{x) indicates a summation over the in-plane nearest 
neighbor lattice sites of x. The quantity b(t;y,y\x) = S yiV '5 s ^b(t;y), which is needed to enforce the condition 

C(t;y,y,0) = 1 for all times t, has to be determined self-consistently^. In addition, the two-time correlator must 
yield the one-time correlator for t = s, i.e. 

C(t,t;y,y',£)=C(t;y,y\x). (17) 

The solution of these equations of motion is outlined in the Appendix. For decorrelated initial conditions, our result 
is: 

d-l 

C(t,s;y,y',x) = (l y - y -(j(t + »)) - I v+y , ( 7 (t + «))) J[ I Xi ( 7 (t + a)) 

i=l 

+ J2 f dTb{T,u)e^ t+s - 2 ^\{l x ^ 1 {t + s-2T)) 

u>l"' i=l 

x (j tt _ tf ( 7 (t - r)) - /u+»(7(* - r))) - t)) - WW* - r))) (18) 

where the functions /„ are modified Bessel functions^ and where we have taken into account the special form of 
b(t;y,y' ,x) and Equation (|A18[) . It remains to fix the parameter b(t,y), which we determine from the condition 
C(t; y, y; 0) = 1. For large d the factor 7=1/ (Td) becomes small, and we can use the following approximation 

W(7(«-t))*V*' + 0(T) (19) 



and similarly for other terms, see also Ref^l. This yields for vanishing layer magnetization the equation 

t 

d T e- 2 ^-^b(r,y) (20) 



for all y and t. This equation can be solved by Laplace transform, yielding the result b(t, y) = 2 for all y and t. It 
then follows that the correlation function in the semi-infinite model is given by Eq. (fTS)) with b(t, u) set to 2. One 
can get rid of the sum over u by using J2m=-oo Im+k{zi)Im{zn) = h{zi + z 2 ) and I n (z) = 7_ n (z). After doing so, we 
obtain for the surface autocorrelation function the expression 

d-l 



C 1 (t,s)=e-^(l ( 1 (t + s))) ~ (j ( 7 ( t + s ))-/ a ( 7 (t + a ))) 
+2 f dre-( f+s - 2T ) (l ( 7 (i + 8 - 2r))) (l Q ( 7 (t + a - r)) - J 2 ( 7 (i + a - 2r))) . (21) 



It is worth noting that in the limit where y sa y' — > 00 we also recover the known bulk behavior of the correlation 
function, as we get with the help of expression (|A18[) the expression 



d-l 



C(t,s;y,y',x) = e-^I y - yl ( 7 (t + a)) ]J / Xi ( 7 (t + a)) 



2 f dr e- {t+s - 2T h y _ yl ( 1 (t + s-2T)) T[l Xi ( 7 (t + s- 2 

Jo i=l 



- 2r)) (22) 
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which is precisely the expression found by Garriga et ali^s in Fourier-space. 

We immediately remark that for a quench inside the ordered phase with T < T c = 1 no simple scaling behav- 
ior of the surface autocorrelation is observed, due to the extremely rapidly (i.e. exponentially) increasing Bessel 
functions. A similar absence of dynamical scaling is also seen in the bulk system quenched below the critical 
pourb2i. At the critical point however, when T = 1 and therefore 7 = 1/d, we can use the approximation 
e~ u I v (u) tn (27TU)" 1 / 2 cxp(— v 2 /{2u)). As the first term in (|2"Tj) decreases more rapidly than the second one, we 
find in the scaling regime (with Y = t/s) 

C 1 (t,s)=4(jf) ^ ' s d ' 2 ({Y -l)- d ' 2 -(Y + ir d / 2 ). (23) 
This allows us to identify both the noncquilibrium exponents b\ and Aci and the scaling function fd(Y) 7 sec Eq. (J5|): 

h = - r \ Cl =d + 2, f Cl (Y) = 4 ((F - l)^ 2 - (Y + l)- rf / 2 ) (24) 

where we used that in the limit of large d the critical dynamical exponent is equal to 2. 



B. The response function 



In order to compute the response function we start from the differential equation 2 ^ 

dt{o2, v {t)) = -{(J S>y {t)) + (t 3 . y {t)) (25) 

for (crg,y(t)) in the presence of a small external magnetic field hf*. As both h e g* and 1/d are small we can develop 
the tanh to first order in both quantities: 

tanhCfcr, V (i)/T) a -hgfrt) + \ (<J S , y+ i(t) + a s , y ^{t) + g ^(t)) . (26) 
where it is understood that c^o = 0. The definition 

ms;y,y';x-x'):= S ^A (27) 

of the response function now directly yields the differential equation (where we set again x ' = 0) 
d t R(t, s; y, y', x) = -R(t, s; y, y', + ~ s; y + 1, y', x) + s; j/ - 1, y', x) 

+ £ R(t, s; y, y', z)) + U(t - s)<W% (28) 



z(x) 



with R(t, s; 0, y' , x) = = R(t, s;y,0, x). This equation is solved with similar methods as outlined in the Appendix 
for the correlation function. As a result we obtain 



d-l 



R(t, s; y, y', x) = ^A e ~^ ( W ( 7 (i s)) - W - *))) LI 7 *« - *)) ( 2 9) 

For the case T = T c = 1 the surface autoresponse function can again be evaluated in the scaling regime, yielding 
(with Y = t/s) 

JfeM = 2 4 ( ^Y^s-i-HY-i)-^ (30) 



T r \ d 

and therefore 



ai = ~ A Rl =d + 2, /fli (y ) = ^^ d/2 ( y (31) 
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We can now also compute the surface fluctuation-dissipation ratio from the expressions (|23|) and (|30[) and obtain 



X!{t,s) 



d a d(t,8) 



(Y-l) 



y((Y - - (Y + l)-*- 1 ) - (V - ~(Y+ 1)-*) 



(32) 



from which the limit value X%° = 1/2 follows. 

Comparing with the results obtained for the spherical model, sec Table HI we note that the values of the nonequilib- 
rium exponents in the critical short-range Ising model in the limit of high dimensions are in full agreement with the 
values obtained for the critical spherical model in dimensions d > Even the scaling functions are identical up to 
a nonuniversal numerical prefactor. Besides, the values of the universal quantities are identical to the values obtained 
in the field-theoretical Gaussian model^. This nicely demonstrates that in the aging regime universal nonequilibrium 
features are indeed encountered close to critical surfaces. 



III. QUENCHING SEMI-INFINITE SYSTEMS FROM HIGH TEMPERATURES: NUMERICAL 

RESULTS 



In the following we present the results of extensive numerical simulations of the standard two-dimensional semi- 
infinite Ising model with only nearest-neighbor interactions quenched inside the ordered phase. In the bulk case this 
model is known to render rather faithfully the physics of real systems undergoing phase ordering. The Hamiltonian 
is given by the usual expression 

H = — CiCj , (33) 

(id) 

where i and j label the sites of a semi- infinite lattice. The sum extends over nearest neighbor pairs, and we have the 
same coupling strength J > for every bond connecting neighboring spins. This system exhibits a continuous phase 
transition at the bulk critical point T c = 2/ln(v / 2 + 1) ps 2.269 (where the temperature is measured in units of J/fc_e, 
with ks being the Boltzmann constant). 

Whereas surface aging behavior has already been studied in the past for Ising models quenched onto the critical 
pourfcH, this does not seem to be the case for quenches below the critical point. The following numerical study 
therefore allows us to close a gap in our understanding of the nonequilibrium dynamical behavior of classical spin 
models. Especially, it yields new insights into the local dynamical behavior of systems undergoing phase-ordering in 
the presence of surfaces. 

For these simulations we use periodic boundary conditions in one direction and free boundary conditions in the 
other direction. We thereby consider square systems with N = L x L spins were L ranges from L = 300 to L = 1000, 
thus making sure that the data obtained at any one of the two surfaces are representative of the semi-infinite system. 
Only data free of finite-size effects are discussed in the following. Our focus lies on the surface autocorrelation function 
and on the surface autoresponse function. The surface autocorrelation function is given by the expression 

C 1 {t,s) = ± (*•(*)<*(«)>, ( 34 ) 

i£surface 

where the sum is over all the spins in the two surfaces. The data discussed in the following have been obtained after 
averaging over at least 5000 different runs with different realizations of the noise. In order to study the response to 
a magnetic field, we apply a weak binary random field between the time t — (at which the quench takes place) 
and the time t = After the field has been switched off, we monitor the decay of the surface thermoremanent 
magnetization given by the expression 

Ml(M) = 2T ^ Jh^M/T (35) 

i^surface 

where hi is the strength of the binary random field at site i. In addition to averaging over the realizations of 
the noise we also average over the realizations of the random field as indicated by the bar. We discuss here data 
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obtained with \hi\ = 0.1 (wc checked that our conclusions remain the same when we slightly vary the value of \hi\). As 
response functions are very noisy, we average over many more runs than for the autocorrelation. The thermoremanent 
magnetization data discussed in this section have been obtained after averaging over typically 200,000 runs. 



A. Autocorrelation function 

Before discussing the surface autocorrelation function, let us briefly mention some results obtained for the autocor- 
relation function in the corresponding two-dimensional bulk system. The expected scaling form 

C(t,s) = s- b f c (t/s) with f c (t/s) ~ (t/s)- Xc/z for t/s > 1 (36) 

has been verified in various numerical studies. These studies showed that 6 = and yielded the value Ac / z = 
0.63(1) & 34 ' 35 (recall that z = 2) for the exponent governing the long-time decay of the scaling function. Numerous 
theoretical approaches have been proposed for computing the scaling function /( g 27 i 28 ' 29 i 30 ' 31 , the most successful 
being the recent exploitation of space-time symmetries within the theory of local scale invarianc o 10 ! 32 ' 33 . 

The main question we address here concerns the scaling behavior of the surface autocorrelation function. Let us 
start by looking at the long-time decay of C\(t, s) with s = 0, as it is well known that this quantity is usually the most 
appropriate for the determination of \c t ■ In Figure 1 we show this quantity for two different temperatures, T = 1 
and T = 1.5, lower than the critical temperature. For comparison we also include the bulk autocorrelation function 
C(t, s = 0) for the same two temperatures. Whereas at short times the surface autocorrelation (this is also true for 
the bulk quantity) is clearly temperature dependent, at longer times the two curves get identical. Interestingly, the 
decay of the surface correlations follow a power-law at late times. This power-law decay is faster at the surface than 
inside the bulk, yielding the value Xci/ Z = 0.95(3) which should be compared to the value \c/ z = 0.63(1) obtained 
inside the bulk. Obviously, this faster decay is due to the reduced coordination number at the surface. 

In Figure 2 we discuss the behavior of the surface autocorrelation function C\(t,s) with s > 0. When plotting 
C\(t, s) versus t/s, we do not observe a data collapse, see Figure 2a, in contrast to the data collapse observed when 
plotting the bulk autocorrelation as a function of t/s. The data shown in Figure 2a at first look suggest that the 
local exponent b\ is different from zero at the surface. A more thorough analysis reveals however that a good scaling 
behavior can not be achieved with a constant b\ > 0. Figure 2b shows our best result obtained for b\ = 0.13. A 
reasonable data collapse can be achieved this way for large values of t/s, but scaling breaks down for t/s < 25. Taken 
at face value, this would suggest for the surface autocorrelation function the existence of a large threshold value of 
t/s below which dynamical scaling is not observed. The possible physical mechanism responsible for this threshold is 
far from obvious. A better data collapse can be achieved by allowing the exponent b\ to depend itself on t/s, but a 
non-constant exponent varying as a function of t/s is not supported by any theoretical approach. 

We propose here another interpretation of the numerical data that is based on the recent observation that large 
finite-time corrections can to some extend mask the true scaling behavior of the autocorrelation function in phase- 
ordering systems*. In order to take the existence of finite-time corrections into account, we try to describe our data 
by the ansatz 

CiM = f Cl (t/s) + s- b 'g Cl (t/s) , (37) 

where the first term is the expected scaling behavior with b\ = 0, whereas the second term is the finite-time correction 
that is of decreasing importance for increasing values of the waiting time s. This ansatz has recently been used for the 
analysis of the autocorrelation functions in disordered ferromagnets quenched below their critical point^i 3 ^. In Figure 
2c we show C\{ys, s) as a function of s for various values of the ratio y = t/s. The lines show that an excellent fitting of 
the data can be achieved with the extended scaling form (|3T[) with a common value b' = 0.49(1). The scaling function 
fci(t/s), obtained after subtracting off the correction term, is shown in Figure 2d. As the curves for the different 
values of s are not distinguishable on the scale of the Figure, we only show selected points as symbols. The data 
collapse shown in Figure 2d supports our interpretation that the true scaling behavior of the surface autocorrelation 
function is masked by strong finite-time corrections. As a consistency check, we note that the data in Figure 2d 
present for large values of t/s a power-law decay with an exponent 0.95(2), in full agreement with the value of \c\l z 
obtained directly from C±(t, s = 0). Even though our data are perfectly described by Eq. (|37p . we must emphasize 
that we do not yet know why this finite-time correction shows up close to the surface but is not encountered inside 
the bulk. 

Let us end the discussion of the surface autocorrelation function by noticing that the value b' = 0.49(1) of the 
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FIG. 1: (Color online) The surface and bulk autocorrelation functions Ci (t, s = 0) and C(t,s = 0) for two different temperatures. 
At the surface the correlations decay much faster than inside the bulk, yielding the value A^/z = 0.95(3) for the long-time 
power-law exponent, considerably larger than the value Xc/z — 0.63(1) in the bulk case. 



correction term exponent is compatible with 1/2 = \/z. However, we refrain from making the conjecture b' = 1/z 
here without having studied other systems with surfaces (as for example semi-infinite Potts models). 



B. Response function 



Before discussing the surface thermoremanent magnetization M\{t,s), let us again first recall the behavior of the 
corresponding bulk quantity. The bulk thermoremanent magnetization M(t, s) is a temporally integrated response 
function that is related to the response function R(t, s) by the integral 

s 

M(t,s) = J duR(t,u) , (38) 
b 

where the integration is over the whole time interval during which the magnetic field was acting on the system. From 
the scaling form (|4|) of R(t, s), we therefore obtain the scaling behavior 

M(t,s) = s- a f M {t/s) (39) 

for the integrated response. Zippold, Kiihn, and Horner— were the first to point out the existence of a sublcading 
correction term which can be quite sizeable. For the thermoremanent magnetization this leads to the following more 
complete scaling behavior—, 



M(t, s) = s~ a f M (t/s) + s- x */ z g M (t/ S ) 



(40) 
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FIG. 2: (Color online) Discussion of the surface autocorrelation Ci(t, s) obtained after quenching the semi-infinite two- 
dimensional Ising model to T = 1. (a) Autocorrelation as a function of t/s for different values of the waiting time s (the 
lowest curve shows the data obtained for the largest value of s). The expected data collapse with bi = is not observed, (b) 
Plotting s°- ri Ci(t,s) versus t/s leads to a collapse of data for large values of t/s, but no scaling is observed for smaller values 
of t/s. (c) Plot of Ci(ys, s) as a function of s for various values of y = t/s. The full lines are fits to the extended scaling form 
(|37l) with b 1 — 0.49. (d) Scaling function /ci (V s ) obtained from the data shown in (a) after subtracting off the finite-time 
correction term. In (c) and (d) error bars are smaller than the symbol sizes. 



The second term in this equation is in fact the response of the system to fluctuations in the initial state, where 
the scaling function gM(t/s) is expected to be proportional to the power-law (t/ s)~ Xr / z — . For the two-dimensional 
Ising model we have a — 1/z = 1/2 and \r/z = 0.63. Therefore this correction to scaling can not be neglected 
but must be included in order to obtain the correct description of the scaling behavior of the bulk thermoremanent 
magnetization*^. 

In Figure 3 we summarize our findings for the surface thermoremanent magnetization in the two-dimensional semi- 
infinite Ising model quenched below the critical point. Figure 3a shows the behavior of this local response as a 
function of t/s for various values of the waiting time s. In a first attempt, we might try to achieve a scaling behavior 
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by assuming that 

Mi(M)=«-°7M 1 (i/s), (41) 

thereby neglecting any possible corrections to scaling. A reasonable scaling behavior is achieved this way for a value 
of ai ~ 0.40, slightly lower than the expected value 1/z = 1/2. For a more thorough analysis we can fix y = t/s and 
plot the response as a function of the waiting time in a log- log-plot. Fitting a straight line to the data, we obtain 
from the slope of that line a value of ai for every considered value of t/s. Thus, we obtain a\ = 0.38 for t/s = 5, 
ai = 0.39 for t/s = 10, a\ = 0.40 for t/s = 15, and a\ = 0.42 for t/s = 20. This points to the existence of a correction 
term that vanishes for increasing values of t/s. In Figure 3b we test the more complete scaling form 

Mi (t,s) = s~ ai f Ml (t/s) + s- x "i / z g Ml (t/s) (42) 

where the correction term with the scaling function gM 1 (t/s) = r\(t/ s)~ Xr i/ z describes the response of the surface 
to fluctuations in the initial state. Plugging in the value X^/z = 0.95 (where we assume that = Ad holds), we 
obtain a consistent description for any t/s with common values r\ = —0.106(1) for the amplitude of the correction 
term and a\ = 0.50(1) for the exponent of the leading term. The correction term being now completely fixed, we 
can subtract it off from the numerical data and obtain the data collapse shown in Figure 3c. Thus, as for the 
thcrmorcmancnt magnetization in the bulk^Jil, we are able to identify the leading correction term and in addition 
obtain the value a\ = a = 1/z. 

We close this section by a brief discussion of the surface fluctuation-dissipation ratio. In Figure 3d we plot the ratio 

TM,(t, s) 

Zl{t ' s) = Tc^s) (43) 

as a function of s/t for various values of s. This ratio yields asymptotically the limit value Af° of the fluctuation- 
dissipation ratio ([H]), as 

Af° = lim ( lim Z^s)) . (44) 

s >oo \t >oo / 

For a given value of the waiting time, the ratio Z\ (t, s) converges towards a constant finite value when s/t — > 0. At first 
look this might seem surprising as in coarsening systems one expects the limit value = 0. However, this constant 
decreases for increasing values of s. Taking into consideration the leading scaling behaviors of C\(t 1 s) ~ fc t (t/s) and 
of Mi(t, s) ~ s" 1 / 2 /mi(^/ s) found in our study as well as the fact that the scaling functions fc^t/s) and /mi(^/s) 
display for large arguments a power-law behavior with the same exponent 0.95, we find that the saturation value 

lim f >oc Z\(t, s) should vanish as s" 1 / 2 . This is indeed verified in Figure 3e, where s x l 2 Z\(t, s) leads to a collapse of 

the data onto a common curve for s/t small. This is also an a posteriori check that we have indeed correctly identified 
the leading scaling behaviors of both the surface autocorrelation and the surface integrated response functions. 



IV. CONCLUSIONS 



In this paper we have extended the investigation of surface aging phenomena to cases not studied previously. On the 
one hand we have computed nonequilibrium surface quantities in the exactly solvable short-range semi-infinite Ising 
model in the limit of a large number of space dimensions, on the other hand we have presented numerical simulations 
of the standard semi-infinite Ising model quenched inside the ordered phase. 

For a quench to the critical point, we added the semi- infinite Ising model in high dimensions to the list of exactly 
solved models. The universal nonequilibrium surface quantities obtained in this study agree with those obtained for 
the critical semi-infinite spherical model^, as expected for a mean-field like model. In Table I we summarize the 
known results for surface aging phenomena in critical systems at the ordinary transition (the only situation studied in 
this paper) by listing the values of the different universal nonequilibrium exponents as well as those of the asymptotic 
value of the fluctuation-dissipation ratio. It is worth mentioning that the existing numerical data for the semi-infinite 
Ising modet^ indicate a non-monotonic behavior of the limit value of the surface fluctuation-dissipation ratio as a 
function of the dimensionality of the system (being ^ for d > 4, then increasing to 0.59 in three dimensions, before 
decreasing to 0.31 in the two-dimensional system). This behavior is unexpected, and a satisfactory explanation is still 
lacking. 
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FIG. 3: (Color online) Discussion of the surface thermoremanent magnetization M\(t,s) obtained after quenching the semi- 
infinite two-dimensional Ising model to T = 1 where a random magnetic field of strength h = 0.1 is applied between t = and 
t — s. (a) Mi(t, s) plotted against t/s for various waiting times s (the lowest curve shows the data obtained for the largest 
value of s). (b) Plot of Mi(ys, s) as a function of s for various values of y — t/s. The full lines are fits to the extended scaling 
form (|42[l . A consistent description of the data for any value of t/s is achieved for Ar/z = 0.95, a = 0.5, and ri = —0.106, see 
main text, (c) Scaling function fM 1 {t/s) obtained from the data shown in (a) after subtracting off the finite-time correction 
term. Error bars are smaller than the symbol sizes, (d) Plot of Zi(t,s), see Equation (|43|l . versus s/t for different waiting 
times, (e) Plotting s 1 ^ 2 Zi(t, s) versus s/t leads to a data collapse for small values of s/t. 



We also presented large-scale numerical simulations of the two-dimensional semi-infinite Ising model undergoing 
coarsening. From these results we conclude that the result b± ^ found for the spherical model^ is not generic but 
that it is very probably an artifact of that rather artificial model. Indeed, the numerical simulations of the more 
realistic two-dimensional Ising model yield b\ = 0. This indicates that generically the exponent b\, that governs the 
scaling of the surface correlations, vanishes, similarly to what is observed inside the bulk. 

One of the main conclusions of our work is that surface aging phenomena in systems undergoing phase-ordering 
display the same general features as bulk aging phenomena. Simple scaling forms prevail asymptotically for two- 
time quantities like the surface autoresponse and the surface autocorrelation functions, and universal nonequilibrium 
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TABLE I: Available values of nonequilibrium critical surface quantities at the ordinary transition determined in aging systems 
quenched to the critical point. 



model 


oi = bi 


^Rx — Ad 




spherical model (2 < d < 4) — 


d 
2 


3d 
2 


1-1 


spherical model (d > 4}^ 


d 
2 


d + 2 


i 

2 


Gaussian model^ 


d 
2 


d + 2 


1 
2 


Ising model in large dimensions 


d 
2 


d + 2 


1 
2 


Ising model in d = 3, ordinary transition^ 


1.24(1) 


2.10(1) 


0.59(2) 


Ising model in d — 3^ 


0.46(1) 


1.09(1) 


0.31(1) 



TABLE II: Available values of nonequilibrium surface quantities determined in aging systems quenched below the critical point. 



model 


ai 


6i 


A«! = Aci 


spherical modeU^ 


d 
2 


1 


| + 2 


Ising model in d = 2 


1 
2 





1.90(6) 



quantities, with values that differ from the bulk values, can also be identified in semi-infinite coarsening systems, 
see Table 2. For finite times, corrections to scaling can be rather important and might even mask the leading 
scaling behavior. In our study of the two-dimensional semi-infinite Ising model we not only identified a sub-leading 
contribution to the thermoremanent surface magnetization (a similar correction also appears inside the bulk) , but we 
also showed the existence of corrections to scaling in the surface autocorrelation function. The physical origin of this 
last term is not yet clear. It is however worth noting that a similar correction term has recently been shown to exist 
for the random bond Ising model quenched below the critical points. 

The semi-infinite geometry discussed in this paper is of course only a special case of a more general wedge-shaped 
geometry. Wedges in critical systems have been studied quite intensively in the pas t 18 i 39 ' 40 i 41 , as they lead to static 
critical quantities whose values depend on the opening angle of the wedge. However, the local critical dynamical 
behavior in a wedge-shaped geometry has not yet been discussed in the literature. Phase-ordering in wedges can 
also be viewed as being one of the simplest cases of phase-ordering in confined geometries. The study of edge aging 
phenomena is therefore the next logical step in the study of local nonequilibrium dynamical behavior in confined 
geometries, and work along this line is in progress. 
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Appendix A. Ising model in high dimensions: Computation of the correlation function 

In this Appendix we compute the one-time and two-times correlation functions for the semi-infinite Ising model in 
high dimensions. We thereby start by defining the operator rjj^ ' , - for a function / : R X Z>o X Z>o x 7L d ~ x — > R: 

n S,s,',*/(*; V, V', x) := d t f(t; y, y', x) + 2f(t; y, y', x) - ^ (/(i; y + 1, y', x) + f(t; y - 1, y', x) 

+ f(t; y, y' + 1,3) + f(t; y,y' - 1, x) + 2 ^ f(t; y, y', f')) , (Al) 

x' (x) 
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where x'(x) denotes the nearest neighbors of x in the layer y. With this, Equation IT4|) reads 

n( t]y,y>, s c(t;y,y',Z) = K*v,i/,Z)- ( A2 ) 

In order to solve this equation we look for the Green's function satisfying the equation 

u tw,y',sd{t\ u , V, v, y', x) = 5(t)5 y ^6 y , tV 5 Sfl (A3) 

and the boundary conditions 

g(t; u, 0, v, y\ x ) = = g(t; u, y, v, 0, x). (A4) 

We can solve this equation by using a Fourier-Sine transformation on y and y' (which makes sure the boundary 
conditions hold) and a normal Fourier transformation on x: 

g(t;u,k,v,k',q) = 2J X! sin ( fc V) sin ( fc ' y') elX u , V, v, y' , x) , (A5) 

y,y'>0 x 

g(t;u,y,v,y',x) = [ f ^- [ d f sin(fc y) sm(k'y') e - lS ^g(t; u, k, v, k', q) (A6) 

Jo n / 2 Jo J B (2n) a 

Here, the sums are over all lattice sites, whereas B = [— ir, tt}' 1 " 1 is the first Brillouin zone and dq = Yit=i dq%- It is 
straightforward to work out equation (|A3j) in Fourier space, which yields 



dtf}(t; u, k, v, k , q) + ui(k, k , q) g(t; u, k, v, k , q; t) = sin(ufc) sm(vk')6(t) (A7) 
where the expression u>(k,k' ,q) is given by 

co(k,k',q) =u(k)+u(k')+u(q) (A8) 

with 

d-l 



w(k) = Q - 7cos(fc)^ , w(fc') = Q -7cos(fc')J , u>(q) =^ Q -27cos(& 



(A9) 



Equation (|A7|) is readily solved and yields the result 

g(t; u, k, V, k' , r; t) = Q(t) s'm(u k) s'm(v k') exp (— u>(k, k', q)t) (A10) 



where Q(t) is the Hcaviside step function. This expression still has to be brought back to direct space using (|A6 
With the integral d k exp(irfe + cos(fc)z) = 2irl r (z) this yields the result 



d-i 



g(t;u,y,v,y',£) = e(t) e - 2 *(/„_ a ( 7 t) - I u+y (~/t))(I v - y ,( 7 t) - W(7*)) \[ I Xi (2~/t). (All) 

z=l 

With the help of this function, the inhomogencous differential equation (TTJJ is solved by 

/>oo 

C(t;y,y',x) =C h (t;y,y',x) + V V / drg(t - r;u,y,v,y' ,x - x')b(T;u,v,x") (A12) 

.. „.^n SI JO 



u,v>0 x' 

where Ch(t;y,y',x) is an arbitrary solution of the homogeneous equation Dj^ , ~Ch(t;y,y',x) = that satisfies the 
boundary condition (|16[) . 

In order to obtain this special solution Ch{t\ y, y' , x) we start from equation (|A3|) with a vanishing right hand side, 

n t .y t y>^Ch{t] y, y', x ) = (A13) 

for which we obtain the solution 

C h (t;k,k',q) = exp(-uj(k,k',q)t)C h (0;k,k',q) (A14) 
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in Fourier space where the initial value C/j(0; k, k! , q) is given by 

4(0; k, k', $ = E E sin ( M fc ) sin ( w fc ')e" i5 '* C(0; u, «, £ ')• (A15) 

u,t>>0 a;' (if) 

The equation (|T5|) for the two-time correlator is also solved in a similar way. We go to Fourier space and get the 
solution 

C(t, s; k, k', q) = exp (-1 (w(jfc) + uj(q)) (t - s)) C(s; k, k',q). (A16) 

Inserting the Fourier transform of the solution of equation (IA12[) . we obtain after transforming back to real space the 
final result 

d-i 

C(t,s;y,y',x) = £ ]T C(0; u, v, x ') e-( t+s > [] I Ti . r[ ( 7 (t + a)) 

x (/ u _ 2/ (7t)-/„+ 3/ (7*))(W(7s)- WM) (A17) 

(/u- y (7(* - r)) - J«+»(7(« - t))) (j«-»(7(» - I")) - Iu+y{l{s - r))) . 

For decorrelated initial conditions C(0;ii, v, x ') = S U _ V S S , g, one can rearrange the Bessclfunctions in the first sum 
using J27=-oo Iu{zi)Iu+k{z2) = h(zi + z 2 ) which gives 

E - *u+»(7*)) (W(7«) - W(7*)) = W(7(< + »)) - W(7(* + «)). (A18) 

«>o 

an equation needed for deriving the final result (|2ip . 
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